CZ_COVID-19 Projekce umrti
Philip Ridgill
10/01/2020
Long-term scenarios of daily deaths
CZ_expected_deaths_timeline_projected_end <- CZ_expected_deaths_timeline_projected %>%
filter(date_expected_death == date_view_model_end)
CZ_expected_deaths_timeline_projected_summary <- CZ_expected_deaths_timeline_projected %>%
mutate(month = month(date_expected_death),
year = year(date_expected_death)) %>%
group_by(year, month, Scenario) %>%
summarise(expected_deaths = sum(expected_deaths)) %>%
mutate(date = paste("15/", month, "/", year),
date = dmy(date))
ggplotly(
CZ_expected_deaths_timeline_projected %>%
left_join(CZ_all, by = c("date_expected_death" = "date")) %>%
filter(date_expected_death <= date_view_model_end) %>%
filter(date_expected_death >= date_view_model_start) %>%
ggplot(aes(x = date_expected_death, y = round(expected_deaths_7,2), col = Scenario)) +
geom_line(size = 1) +
geom_text_repel(aes(label = round(expected_deaths_7,0), colour = Scenario), data = CZ_expected_deaths_timeline_projected_end, size = 3, vjust = 0, hjust = -200) +
#geom_line(aes(y = expected_deaths), col = "grey10") +
geom_line(aes(y = new_deaths_7), col = "black", size = 1) +
geom_vline(xintercept = date_model, col = "red", linetype = "dashed") +
scale_color_manual(values = c("#4DAF4A", "#377EB8", "purple", "#E41A1C")) +
scale_x_date(date_breaks = "1 month", date_minor_breaks = "1 week", date_labels="%b") +
theme_light() +
theme(panel.grid.minor = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend(nrow=2,byrow=TRUE)) +
labs(y = "Deaths (7 day rolling average)", x = "", title = "Scenarios of expected daily deaths", subtitle = "Based on past data and 4 scenarios of growth in cases")
) %>%
layout(legend = list(
orientation = "h",
x = -0,
y = -0.1
)
)CZ_expected_deaths_timeline_projected_end <- CZ_expected_deaths_timeline_projected %>%
filter(date_expected_death == date_view_model_end)
CZ_expected_deaths_timeline_projected_summary <- CZ_expected_deaths_timeline_projected %>%
mutate(month = month(date_expected_death),
year = year(date_expected_death)) %>%
group_by(year, month, Scenario) %>%
summarise(expected_deaths = sum(expected_deaths)) %>%
mutate(date = paste("15/", month, "/", year),
date = dmy(date))
CZ_expected_deaths_timeline_projected %>%
left_join(CZ_all, by = c("date_expected_death" = "date")) %>%
filter(date_expected_death <= date_view_model_end) %>%
filter(date_expected_death >= date_view_model_start) %>%
ggplot(aes(x = date_expected_death, y = round(expected_deaths_7,2), col = Scenario)) +
geom_line(size = 1) +
geom_text_repel(aes(label = round(expected_deaths_7,0), colour = Scenario), data = CZ_expected_deaths_timeline_projected_end, size = 3, vjust = 0, hjust = -1) +
# geom_text(aes(x = date, y = 0, label = round(expected_deaths,0), colour = Scenario), data = CZ_expected_deaths_timeline_projected_summary %>% filter(Scenario == "My projection") %>% filter(date <= date_view_model_end) %>%
# filter(date >= date_view_model_start), size = 3, vjust = -1, hjust = 0) +
#geom_line(aes(y = expected_deaths), col = "grey10") +
geom_line(aes(y = new_deaths_7), col = "black", size = 1) +
geom_vline(xintercept = date_model, col = "red", linetype = "dashed") +
scale_color_manual(values = c("#4DAF4A", "#377EB8", "purple", "#E41A1C")) +
scale_x_date(date_breaks = "1 month", date_minor_breaks = "1 week", date_labels="%b") +
theme_light() +
theme(panel.grid.minor = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend(nrow=2,byrow=TRUE)) +
labs(y = "Deaths (7 day rolling average)", x = "", title = "Scenarios of expected daily deaths", subtitle = "Based on past data and 4 scenarios of growth in cases")
Model description
1) Divide population into age cohorts
plot_cases <- CZ_positive_dead_age_category %>%
ggplot(aes(x = age_category, y = count_cases)) +
geom_col(fill = "red") +
theme_light() +
theme(panel.grid.minor = element_blank(),
axis.title.x = element_blank(),
legend.position = "none")
plot_deaths <- CZ_positive_dead_age_category %>%
ggplot(aes(x = age_category, y = count_dead)) +
geom_col(fill = "black") +
theme_light() +
theme(panel.grid.minor = element_blank(),
legend.position = "bottom") +
labs(x = "Age category")
grid.arrange(plot_cases, plot_deaths, nrow = 2, ncol=1)
2) Calculate fatality per age group
Delay the date of positive cases by 6 to 14 days (depends on age) to reflect the delay between a positive test and death and get a better estimate of fatality (deaths/cases)
(plot_fatality <- CZ_fatality %>%
ggplot(aes(x = age_category, y = fatality, label = paste(100*round(fatality,4),"%"), fill = -fatality)) +
geom_col() +
geom_text(vjust = -0.4) +
scale_y_continuous(labels = percent_format()) +
theme_light() +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
legend.position = "none")
)
3) Use an estimate of the delay between a positive test and death
(plot_death_timing <- CZ_death_timing %>%
ggplot(aes(x = delay, y = probability)) +
geom_col() +
theme_light() +
scale_x_continuous(breaks = seq(0,25,2)) +
facet_wrap(~ age_category, ncol = 2) +
theme(panel.grid.minor = element_blank()) +
labs(title = "Timing of death by age category", y = "Probability", x = "Number of days since positive test")
)
Backtest: Check estimate looks accurate
Expected deaths = number of cases * fatality * probability of death occuring that specific day
CZ_expected_deaths_timeline_by_age %>%
left_join(CZ_deaths_by_age, by = c("date", "age_category")) %>%
filter(date <= date_today - 7) %>%
ggplot(aes(x = date, y = expected_deaths_7, group = age_category)) +
geom_line(col = "brown") +
#geom_line(aes(y = expected_deaths), col = "grey10") +
geom_line(aes(y = deaths_7), col = "black") +
geom_vline(xintercept = date_today - 7, col = "red", linetype = "dashed") +
theme_light() +
facet_wrap(~ age_category) +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank()) +
labs(y = "deaths (7 day rolling average)", x = "date", title = "Expected deaths (brown) vs actual deaths (black)", subtitle = "Cut off 7 days ago due to reporting lags")
4) Create scenarios of future growth in cases
Use week on week growth of cases as a simple R measure. These R values are entirely arbitrary.
CZ_cases_by_age_projected %>%
filter(date >= as.Date("01/01/2021", format = "%d/%m/%Y")) %>%
select(week, Scenario, R) %>%
filter(week >= week_now,
!is.na(R)) %>%
distinct() %>%
mutate(R = R^7) %>%
ggplot(aes(x = week, y = R, group = Scenario, col = Scenario)) +
geom_line(size = 1) +
scale_color_manual(values = c("#4DAF4A", "#377EB8", "purple", "#E41A1C")) +
scale_x_continuous(breaks = seq(1, 53, 1)) +
theme_light() +
theme(panel.grid.minor = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend(nrow=2,byrow=TRUE)) +
labs(y = "Simplified R", x = "Weeks", title = "R of different scenarios")
5) Calculate expected deaths for each age group
For each positive case apply the respective fatality for their age group and multiply by the probability of the person dying each day. For example if a 90+ year old has a 20% fatality rate and I estimate 20% of those who die will die after 6 days, the overall probability of that person dying (expected death) on the date of infection + 6 will be 0.2*0.2 = 0.04.
Then sum up each person’s probability of dying (expected death) for each day to get the number of expected deaths for each day.
CZ_expected_deaths_timeline_by_age_projected %>%
left_join(CZ_all, by = c("date_expected_death" = "date")) %>%
filter(date_expected_death <= date_view_model_end) %>%
filter(date_expected_death >= date_view_model_start) %>%
filter(Scenario == "My projection") %>%
ggplot(aes(x = date_expected_death, y = round(expected_deaths_7,2), fill = age_category)) +
geom_col() +
geom_line(aes(y = new_deaths_7), col = "black", size = 1) +
geom_vline(xintercept = date_model, col = "red", linetype = "dashed") +
scale_x_date(date_breaks = "1 month", date_minor_breaks = "1 week", date_labels="%b") +
theme_light() +
theme(panel.grid.minor = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend(nrow=2,byrow=TRUE)) +
labs(y = "Deaths (7 day rolling average)", x = "", title = "Expected daily deaths by age category", subtitle = "Based on past data and my main scenario of growth in cases")
6) Finally sum up each age category to get the overall projection
The End! Sort of, review output corresponds to past reality -> if not adjust model.
7) Additional step for “complex” model
To adjust for a slightly less accurate projection in the second wave I incorporated % of positive tests (positivity) as a measure of how early/late people get tested. I found that the % of positive tests doesn’t increase fatality (at least according to official covid data), but a high % of positive tests meant people were tested 2-3 days later than usual, meaning people died earlier than my model projected.
To adjust for that I developped a scalar which adjusts the probabilities of the delay between people being tested and dying. In simpler terms if positivity is > 20%, people die earlier than in the baseline model and if positivity < 20% people die later. The scalar hits its maximum level (adjustment) when positivity hits 33% -> adjusts delay by 2-3 days.
However, recently the introduction of AG testing seems to alter the relation between positivity and delay to testing as people are being tested “on time” despite a higher PCR positivity rate.